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Abstract 



We introduce an approach for performing "very long" computer 
simulations of the dynamics of simplified, folded proteins. Using an 
alpha-carbon protein model and a fine grid to mimic continuum com- 
putations at increased speed, we perform unbiased simulations which 
exhibit many large-scale conformational transitions at low cost. In the 
case of the 72-residue N-terminal domain of calmodulin, the approach 
yields structural transitions between the calcium-free and calcium- 
bound structures at a rate of roughly one per day on a single Intel pro- 
cessor. Stable intermediates can be clearly characterized. The model 
employs Go-like interactions to stabilize two (or more) experimentally- 
determined structures. The approach is trivially parallelizable and 
readily generalizes to more complex potentials at minimal cost. 

1 Introduction 

The biological functions of many proteins result from their folded-state dy- 
namics, and in particular, from conformational changes among meta-stable 
states. Motor proteins perform their vital functions based upon structural 
transitions, and additionally, such transitions often accompany ligand-binding 
and catalysis events: those occurring in calmodulin and in adenylate ki- 
nase upon ligand binding are textbook examples pQ among numerous others. 
However, little is known about the detailed dynamics of these large-scale 
transitions. How cooperative are the changes? Do meta-stable intermediate 
states often occur? How many reaction pathways play an important role? 
Do any generic kinetic or dynamical features appear, especially in light of 
the expectation that nature tunes barrier heights to a few times the thermal 
energy scale ksT ~ RT1 

Calmodulin (CaM), a 148-residue calcium-binding and signalling protein, 
is an ideal test system because of its modest size and because it has been 
the subject of intensive experimental scrutiny (e.g., (3131111301111313 
ITU]). CaM undergoes large scale conformational transitions both within and 
between its two domains. CaM has also been studied computationally, for 
instance by the Garcia and Kuczera groups [TTJ H3 [T3] , but these all-atom 
molecular dynamics studies have been limited to nsec timescales. 

Many approaches have been developed for the general problem of deter- 
mining reaction pathways — including both "static," "quasi-dynamic," and 
"ensemble" approaches. We term "static" those approaches which yield a 
single, presumably optimal pathway (e.g., [HI HSJ ifel HZl UHl UHl EOl EH E2J 
E3 [24J ) . The so-called static picture, by definition, does not account for 
inherent thermal fluctuations or the possibility of multiple pathways, though 
Elber and Shalloway modelled the effect of temperature in a static picture 
|25| . Biased molecular dynamics approaches, like targeted and steered dy- 
namics [23 123 123 I2H| attempt to include more realistic aspects of the true 
dynamics; such methods typically generate a small number of potentially 
dominant dynamical pathways in protein systems, rather than a full ensem- 
ble. 



Figure 1: A schematic bead model of a protein on a fine grid, and the residue- 
level depiction of calmodulin's N-terminal domain. The centers of the residue 
beads (alpha-carbon locations) are restricted to grid sites. Because the grid 
spacing is much less than the size of the beads, the discretized simulation 
mimics a continuum calculation. Interaction potentials for all possible bead 
orientations are stored in memory during the simulation, saving potentially 
costly calculation. The right panel depicts an averaged simulation structure, 
starting from the configuration of PDB code lcfd. 



"Ensemble" approaches attempt to generate suitably distributed sets of 
reacting trajectories. The pioneering work of Pratt employing the Metropolis 
approach [3U] was applied by Chandler and coworkers to a variety of problems 
(e.g., |3TJ E21 ESI EH EI3 EE| ) • Elber and coworkers developed the stochastic 
difference equation approach, which generates approximate reactive trajec- 
tories without high-frequency motions [37J EH EH]- Woolf and Zuckerman 
purused a non- Metropolis ensemble approach [1UJ |HJ E21 US] , as did Mazonka 
et al. jH] ; see also the work of Eastman et al. [IS] . A notable, independent 
approach for studying rare events is the "weighted-ensemble Brownian dy- 
namics" of Huber and Kim j^B] , which is conceptually related to the earlier 
work of Harvey and Gabb [26j. 

One recent study pursued similar goals to the present investigation for 
a cubic-lattice toy model. Specifically, Borovinskiy and Grosberg studied 
the design of cubic-lattice "protein" models capable of undergoing conforma- 
tional transitions |T7j . 

Here we introduce an unbiased methodology to study large-scale struc- 
tural transitions in folded proteins, and we use it to perform an initial exam- 
ination of the structural dynamics in the N-terminal domain of calmodulin. 
The method builds on two existing computational strategies. The first, due 
to Panagiotopoulos and Kumar jJHIESj, is a general discretization approach 
which involves the use of a fine grid (with lattice spacing much less than the 
particle or atom size) to mimic continuum calculations. Atoms or particles 
are only permitted to occupy grid sites, and interaction energies or forces are 



computed in advance for all necessary inter-particle displacements; see Fig. 
^ Panagiotopoulos and Kumar found that critical and coexistence proper- 
ties of fluids could be reproduced very precisely with a speed gain of one to 
two orders of magnitude [JH1IIH]- The notion of using a fine grid to simulate 
biomolecules also conceptually draws from coarse-grained lattice simulations 
of proteins [SUl EU E2J EH1 EH ESI EE1 EI] , particulularly the high coordination 
lattice studies by Kolinski, Skolnick, and coworkers. 

The second principal computational precursor to the present work is the 
Go model [SHIED], which represents proteins in solely a structural way, with- 
out reference to the underlying chemistry. In particular, the Go model em- 
ploys a chain of residue "beads", each of which only reacts favorably with 
other residues which are nearby in the native structure. Non-native contacts 
are penalized. This simple model was developed to study the dynamics of 
protein folding and it is still widely used for that end [SHI EH EI] • Note that 
the Go model is trivially generalized to other levels of chemical detail — e.g., 
to all-atom simulation jo"!] . 

While the present approach is not restricted to simplified models, we 
use a highly reduced (residue-level) model both because it is easier to im- 
plement and because residue-level models have already provided meaningful 
bio-macromolecular results. Reduced models of proteins have often been 
used in the past, particularly to study protein folding dynamics and ther- 
modynamics (e.g., [naE3E3EaEiiEaEniEniEziEHiE3izniiziiiza 

EH D31 E3 E01 113 EE EH), "ab initio" protein strucutre prediction (e.g., 
IZHl|73IHniEIllH2llHSllHllHSlESl), and coarse-grained dynamics (e.g., jH7j); 
see also a review of structure-based potentials [HE]- Moreover, it has been 
shown conclusively that even models which do not distinguish among atom or 
amino acid types can capture intermediate scale (i.e., alpha carbon fluctua- 
tions) and large scale motions (i.e., the slowest modes) [5§] W5[ I§T1 IM] . 

The present, unbiased approach should also prove useful in the refinement 
of the biased "ensemble" methods discusssed above. Because the ensemble 
methods have not been fully vetted in protein systems, the present approach 
can usefully contribute "perfect" ensembles of unbiased large-molecule tran- 
sition trajectories which can serve as "gold standards" for comparison with 
biased methods. In particular, unbiased ensembles can help shift the focus to 
a number of biologically and methodologically important questions regard- 
ing intermediates: (i) can their lifetimes be estimated accurately? (ii) do 
intermediates introduce a very wide variation in the durations of transition 
events? (iii) are multiple pathways typically observed? 

In outline, the next section, [21 describes the simulation approach used 
here. We emphasize that our method is not tied to the choice of model (i.e., 
forcefield), which is highly simplified for this initial study and is described in 
Sec. El Sec. |U presents the simulation results, and highlights the straightfor- 
ward observation of stable intermediates. Our concluding discussion is given 
in Sec. El 
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Figure 2: The highly reduced, Go-like, inter-residue potentials used in the 
simulations. The "non-native" pair potential is used between residues that 
are not in contact in the holo or the apo structures of calmodulin. A contact is 
defined by a distance Rij < 8 A between alpha-carbon atoms of residues i and 
j in the Protein Data Bank coordinate files. The "single-native" potential 
applies for residue pairs in contact in only one of the two structures, while 
"double-native" is for pairs contacting in both. 



2 Methodology 

2.1 Fine-grid simulation 

The fine-grid simulation methodology is adopted directly from Panagiotopou- 
los and Kumar [ISJ H§]- The idea is to allow particles (here, residues) to 
occupy only discrete positions on a grid. Interaction energies are then stored 
in arrays, rather than computed at each dynamics step. 

More specifically, a "fine" grid is used — i.e., one in which the grid spacing 
is much smaller than the particle size — in an attempt to mimic continuum 
results. This point bears emphasizing: the grid approach is adopted solely 
as a means to reproduce continuum results at greater computational speed; 
the errors or artifacts resulting from using a fine grid are expected to be 
negligible due to the small size of the grid spacing. Indeed, the central result 
of Panagiotopoulos and Kumar jlHl HH] is that very sensitive liquid-vapor 
coexistence curves and critical points could be determined quite accurately 
within the grid approach. A rough rule of thumb from Panagiotopoulos' 
work is that a lattice spacing of 1/5 the particle size is sufficient for practical 
computations. 

How fine a grid spacing is necesary in the present context? Owing to 
the crudeness of the models used here (see Sec. EJ), we do not require high 
precision. However, to implement a dynamic Monte Carlo scheme (see below) 
with a sufficient acceptance ratio, a small grid spacing is required. In the 
simulations discussed here, the lattice spacing was 0.13 A, which may be 
compared with residue-residue interaction distances of 4 - 8 A. In other words, 
a very fine grid is employed here. 




Figure 3: The "start" and "end" calmodulin structures used to construct the 
Go-like potential. The left panel is the calcium-free (apo) structure, PDB 
code lcfd, determined by solution NMR, from which the simulations were 
started. The right panel depicts the calcium bound {holo) strucutre, PDB 
code lcll, determined from X-ray diffraction. 



2.2 Dynamic Metropolis Monte Carlo 

In this initial implementation, we have chosen to use dynamic Monte Carlo 
(DMC) dynamics, which is a common, current choice for simulations of pro- 
tein folding dynamics (e.g., [SHI EDI IS])- While DMC dynamics were chosen 
because of the relative ease of implementation, they may be justified on two 
grounds, (i) Because only small, local trial moves will be considered — see 
below — and accepted with a Boltzmann-factor-preserving probability, DMC 
may be considered a variant of overdamped Brownian dynamics. Brownian 
dynamics have a well-understood physical basis [HH]. (ii) Further, due to the 
highly simplified nature of the models used here (see Sec. El), one can expect 
the simulated dynamics — of any kind — only to give a qualitative pic- 
ture. This qualitative description should not be sacrificed by Brownian-like 
dynamics embodied in DMC. 

The dynamic Monte Carlo is designed to make highly local moves in a 
trivial way: at every time step, a randomly chosen residue makes a trial move 
with uniform probability to one of the 26 lattice sites on the surface of the 
3x3x3 cube centered at its present location. The move is accepted or 
rejected according to the usual Metropolis criterion: accept if and only if the 
energy decreases or a random number R u chosen uniformly from the interval 
< R u < 1 is less than exp (-AU/RT), where AU is the change in the total 
energy (Hamilitonian) of the system. 

2.3 Hardware, Software, and Computer Time 

The results reported below are based on a computer program written in 
the C language and compiled with the "gcc" compiler on Linux machines. 



Simulations were run on Intel processors of both 2.4 and 2.8 GHz, which 
permitted approximately 3-4 xlO 8 Monte Carlo steps per hour. All the 
trajectories discussed here were obtained in less than one week of of processor 
time. 

3 Model 

The residue-level model describes a highly reduced protein without any true 
chemistry, following the spirit of many previous workers [2BJ E01 EH HM E3] ■ 
Nevertheless, the model does include the essential toplogical features of 
calmodulin: residue connectivity and sterics, as well as the attractive in- 
teractions found in both apo and holo strucrures. Only the N-terminal do- 
main was included — specifically the 72 residues numbered 4 - 75 in PDB 
structures lcfd (apo, unbound) and lcll (holo, bound to calcium ions). 

The model of the N-terminal domain of calmodulin is specified by cova- 
lent, steric, and attractive interactions; see Fig. El Consectutive residues are 
"covalently" connected by an infinitely deep, attractive square well centered 
on the native-state separation distance of structure lcfd; the well width was 
taken to be 10% of the native separation (±5%). Sterically, each residue ex- 
cludes other residues from its hard-core exclusion zone, a sphere with radius 
set to 47.5% (i.e., 95%/2) of the minimum distance to other non-consecutive 
residues in the lcfd native state. Finally, Go-like attractive interactions are 
represented by square wells (of uniform depth e) centered at residue separa- 
tion distances (±5%) of less than 8 A from either of the native structures (i.e., 
PDB codes lcfd or lcll). All simulations were performed at a temperature 
given by ksT = 0.6e. 

The double-well, attractive interactions depicted in Fig. |21 merit further 
justification. At first glance, such a potential may seem unnatural, particu- 
larly to those readers familiar with the distinction between the pair potential 
energy and potential of mean force (see, e.g., [56.): for a simple fluid consist- 
ing of spherical particles, it is well known that a single-well potential gives 
rise to a multi-well potential of mean force. However, amino acids are not 
spherical, and if one insists on a point-residue model (as here), the only way 
to represent direct interactions of differing orientations is for the potential 
to possess more than one well. A second, perhaps facile, justification comes 
from the results: the double well potential accomplishes the principal aim of 
stabilizing two distinct conformational states. Indeed, it requires roughly a 
day of computer time to generate a fluctuation sufficiently large as to jump 
to the second state. A third justification is more pragmatic: when more 
realistic potentials (e.g., [Ml E23 ESI El E9J ESS EFJ OT are discretized with 
the find-grid approach, presumably the "unnatural" double-Go interactions 
can be reduced to minimal levels required to stabilize the distinct states. 




Figure 4: Two sample structural transitions in the simplified N-terminal 
domain of calmodulin, monitored by the numbers of native contacts in the 
apo and holo states. Contacts are defined by inter-C Q distances of less than 
8 A and native contacts are based on comparison to the PDB structures lcfd 
(apo) and lcll (holo). Each of the two trajectories shown required roughly 
30 minutes of single-processor (2.8 GHz) computer time, and such transition 
events occur roughly once per day on a single processor. 



4 Results 

In this preliminary report, all the results are depicted in Figs. 0] - and 
described in the figure captions. To summarize, sample transition events 
are shown (Fig. EJ), and traces of the helix angles confirm that the transition 
are genuine and sustained (Fig.EJ. Next, long-lived intermediates are readily 
visualized using compound variables which simply measure distances between 
helix ends (Fig. EJ), and finally the dynamical evolution of individual residue 
contacts is examined (Fig. EJ). 
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Figure 5: The two transitions of Fig. |3] confirmed by comparison to the 
inter-helix angles in the PDB structure lcll. Helix axes were defined to be 
the vectors associated with the smallest moment of inertia of the residues 
present in the native structures. 




Figure 6: Long-lived intermediates in the transition events of Fig. The 
central panel shows parametric traces of two independent transition trajec- 
tories, with the axes representing distances between the indicated helix ends; 
note that the N-terminal "start" of the peptide is the red helix. The sur- 
rounding structures are averaged over 10 7 MC steps, but remained localized 
at the indicated positions. Note that the two trajectories take slightly dif- 
ferent pathways in the region where vertical axis ranges between 18 and 
25; the lower-right intermediate appears to be "on-pathway" for the orange 
trajectory but not for the purple. 

5 Summary and Conclusions 

We have introduced an approach for studying the long-time dynamics of 
large-scale conformational transitions in proteins. The technique was tested 
on a highly reduced, residue-level model of the N-terminus of calmodulin, 
which undergoes a dramatic rearrangement of its four helices upon the (un- 
binding of calcium ions. Such rearrangements were observed approximately 
once per day on a single processor in unbiased dynamic Monte Carlo simu- 
lation. The scheme is quite promising because it proved capable of spanning 
the range of timescales from that for inter-residue vibrations to the "waiting 
time" between transition events — i.e., the inverse rate. 

The approach is readily extendable in two important ways. First, it is em- 
barrassingly parallelizable, so an inexpensive, modest-sized Beowulf cluster 
can readily improve upon the single-processor output by an order of magni- 
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Figure 7: Dynamical behavior of residue contacts. The left panel depicts 
the evolution of the contacts of three residues during the transition event of 
trajectory 79 (see Fig. HJ). The ll holo contacts" for residue i is the count of 
all residues which are both within R c = 8 A of residue i at the time point 
shown, and which are also within R c of i in the holo structure, lcfd, depicted 
in Fig. El Smoothed versions of the residue contacts evolution are also shown 
on the left for visual clarity. The right panel shows the average differences 
in holo contacts before and after the transition point at 8.610 7 MC steps, 
for every residue in both transition of Fig. EJ Note the similarity of the 
traces for the two transitions, which — interestingly — are not correlated 
with contact-count differences based on the static PDB structures. 



tude or more. Second, the fine-grid discretization approach |48| H2] — which 
simply mimics continuum calculations — lends itself to more complex po- 
tentials: while the addition of particles (e.g., "beads") to any model entails 
a cost, the use of a more chemically realistic potential will require minimal 
additional overhead. Furthermore, other models — reduced or atomistic - 
can be stabilized in experimentally determined structural states using the 
same Go-like interactins employed here. 

Despite the apparent simplicity of the model employed here, the simu- 
lated structural transitions in calmodulin exhibited highly complex behavior. 
Distinct transition pathways and quasi-stable intermediate states could be 
readily identified, as could critical residues. The transition details identified 
in the present study cannot immediately be identified with expectations for 
experimental outcomes because of the simplicity of the model. Nevertheless, 
the results should capture those aspects of the dynamics governed by ster- 
ics and connectivity, in analogy with the success of network-like models in 
predicting large-scale fluctuations [El EH ED El E3l El ■ 

Future studies with more chemically realistic models (e.g., [SU ESI ESI 
1731 1731 175^ I8TH loTH 1751 [70]) implemented in the rapid, fine-grid approach 
should suggest testable experimental hypotheses. Note that the effects on 
dynamics of varying forcefields have recently been addressed by Freed and 
coworkers jSZj. A related question, which also can be addressed directly in 
the present scheme, regards the effects of dynamics type — e.g., Metropolis, 
Langevin, molecular dynamics — on transition events. 



The future determination of reaction pathways and dynamics in a hierar- 
chy of models will also address a question both scientifically and methodolog- 
icallly critical: What is the minimal model necessary to capture residue-level 
features of structural transitions in proteins ? From the "scientific" perspec- 
tive, a minimal model will suggest which interactions govern the transition 
events, while optimal methods should use that model which is least costly, 
computationally. 

Finally, it is appropriate to ask what the unbiased transitions observed in 
the present study presage for the use of the biased "ensemble" methods dis- 
cussed in the introduction. The observation of multiple pathways and widely 
disparate transition event durations suggests that great care will be required 
in any biased effort to generate ensembles of transition trajectories. Never- 
theless, the present approach should provide critical input for testing biased 
methods in protein systems — namely, unbiased trajectories to be used for 
comparison. From a practical standpoint, a set of reduced- model transition 
events may also provide useful starting points to sidestep the potential for 
trapping associated with the Pratt approach. 

Acknowledgments 

The author has benefitted greatly from discussions with many scientists: Ivet 
Bahar, David Deerfield, Jeffery Evanseck, William Furey, Hagai Meirovitch, 
Robert Swendsen, Dror Tobi. Special thanks are due Professor Bahar for her 
invaluable advice, support, and encouragement. 

References 

[1] J. M. Berg, J. L. Tymoczko, and L. Stryer. Biochemistry (5th ed.). 
Freeman, New York, 2002. 

[2] Y. S. Babu, C. E. Bugg, and W. J. Cook. Structure of calmodulin refined 
at 2.2 angstroms. J. Molec. Bio., 204:191, 1988. 

[3] M. Ikura, G. M. Clore, A. M. Gronenborn, G. Zhu, C. B. Klee, and 
A. Bax. Solution structure of a calmodulin-target peptide complex by 
multidimensional NMR. Science, 256:632-638, 1992. 

[4] R. Chattopadhyaya, W. E. Meador, A. R. Means, and F. A. Quiocho. 
Calmodulin structure refined at 1.7 angstroms. J. Molec. Bio., 228:1177, 
1992. 

[5] W. E. Meador, A. R. Means, and F. A. Quiocho. Modulation of calmod- 
ulin plasticity in molecular recognition on the basis of x-ray structures. 
Science, 262:1718-1721, 1993. 

[6] H. Kuboniwa, N. Tjandra, S. Grzesiek, H. Ren, Klee C.B., and Bax 
A. Solution structure of calcium-free calmodulin. Nature Struct. Bio, 
2:768-776, 1995. 



[7] A. Miyawaki, J. Llopis, R. Heim, J. M. McCaffery, J. A. Adams, 
M. Ikura, and R. Y. Tsien. Fluorescent indicators for Ca-2+ based on 
green fluorescent proteins and calmodulin. Nature, 388:882-887, 1997. 

[8] J. Evenas, S. Forsen, A. Malmendal, and M. Akke. Backbone dynam- 
ics and energetics of a calmodulin domain mutant exchanging between 
closed and open conformations. J. Molec. Bio., 289:603-617, 1999. 

[9] W. S. VanScyoc, B. R. Sorensen, E. Rusinova, W. R. Laws, J. B. A. Ross, 
and M. A. Shea. Calcium binding to calmodulin mutants monitored 
by domain-specific phenylalanine and tyrosine fluoresence. Biophys. J., 
83:2767-2780, 2002. 

[10] A. M. Weljie, A. P. Yamniuk, H. Yoshino, Y. Izumi, and H. J. Vogel. 
Protein conformational changes studied by diffusion NMR spectroscopy: 
Application to helix-loop- helix calcium binding proteins. Protein Sci., 
12:228-236, 2003. 

[11] D. Vigil, S. C. Gallagher, J. Trewhella, and A. E. Garcia. Functional dy- 
namics of the hydrophbic cleft in the N-domain of calmodulin. Biophys. 
J., 80:2082-2092, 2001. 

[12] C. Yang and K. Kuczera. Molecular dynamics simulations of calcium- 
free calmodulin in solution. J. Biomolec. Struc. Dyn., 19:801-819, 2002. 

[13] G. S. Jas and K. Kuczera. Free-energy simulations of the oxidation of 
C-terminal methionines in calmodulin. Proteins, 48:257-268, 2002. 

[14] K. Miiller and L. D. Brown. Location of saddle points and minimum 
energy paths by a constrained simplex optimization procedure. Theoret. 
Chim. Acta (Berlin), 53:75-93, 1979. 

[15] M. Berkowitz, J. D. Morgan, J. A. McCammon, and S. H. Northrup. 
Diffusion-controlled reactions: A variational formula for the optimum 
reaction coordinate. J. Chem. Phys., 79:5563-5565, 1983. 

[16] J. Janin and S. J. Wodak. Reaction pathway for the quaternary structure 
change in hemoglobin. Biopolymers, 24:509-526, 1985. 

[17] R. Czerminski and R. Elber. Reaction path study of conformational 
transitions and helix formation in a tetrapeptide. Proc. Nat. Acad. Sci., 
86:6963-6967, 1989. 

[18] R. Czerminski and R. Elber. Reaction-path study of conformational 
transitions in flexible systems - applications to peptides. J. Chem. Phys., 
92:5580-5601, 1990. 

[19] C. Choi and R. Elber. Reaction path study of helix formation in 
tetrapeptides: Effect of side chains. J. Chem. Phys., 94:751-760, 1991. 



[20] R. E. Gillilan and K. R. Wilson. Shadowing, rare events, and rubber 
bands, a variational Verlet algorithm for molecular dynamics. J. Chem. 
Phys., 97:1757-1772, 1992. 

[21] S. Fischer and M. Karplus. Conjugate peak refinement: an algorithm 
for finding reaction paths and accurate transition states in systems with 
many degrees of freedom. Chem. Phys. Lett., 194:252-261, 1992. 

[22] E. M. Sevick, A. T. Bell, and D. N. Theodorou. A chain of states method 
for investigating infrequent event processes occurring in multistate, mul- 
tidimensional systems. J. Chem. Phys., 98:3196-3212, 1993. 

[23] M. K. Kim, G. S. Chirikjian, and R. L. Jernigan. Elastic models of 
conformational transitions in macromolecules. J. Molec. Graph. Model., 
21:151-160, 2002. 

[24] M. K. Kim, R. L. Jernigan, and G. S. Chirikjian. Efficient generation 
of feasible pathways for protein conformational transitions. Biophys. J., 
83:1620-1630, 2002. 

[25] R. Elber and D. Shalloway. Temperature dependent reaction coordi- 
nates. J. Chem. Phys., 112:5539-5545, 2000. 

[26] S. C. Harvey and H. A. Gabb. Conformational transitions using molecu- 
lar dynamics with minimum biasing. Biopolymers, 33:1167-1172, 1993. 

[27] J. Schlitter, M. Engels, P. Kruger, E. Jacoby, and A. Wollmer. Targeted 
molecular-dynamics simulation of conformational change - Application 
to the T <-> R transition in insulin. Molec. Sim., 10:291-309, 1993. 

[28] J. F. Diaz, B. Wroblowski, J. Schlitter, and Y. Engelborghs. Calcu- 
lation of pathways for the conformational transition between the GTP- 
and GDP-bound states of the Ha-ras-p21 protein: Calculations with ex- 
plicit solvent simulations and comparison with calculations in vacuum. 
Proteins, 28:434-451, 1997. 

[29] S. Izrailev, S. Stepaniants, M. Balsera, Y. Oono, and K. Schulten. Molec- 
ular dynamics study of unbinding of the Avidin-Biotin complex. Bio- 
phys. J., 72:1568-1581, 1997. 

[30] L. R. Pratt. A statistical method for identifying transition states in high 
dimensional problems. J. Chem. Phys., 85:5045-5048, 1986. 

[31] C. Dellago, P. G. Bolhuis, F. S. Csajka, and D. Chandler. Transition 
path sampling and the calculation of rate constants. J. Chem. Phys., 
108:1964-1977, 1998. 

[32] C. Dellago, P. G. Bolhuis, and D. Chandler. Efficient transition path 
sampling: Application to Lennard- Jones cluster rearrangements. J. 
Chem. Phys., 108:9236-9245, 1998. 



[33] F. S. Csajka and D. Chandler. Transition pathways in a many-body 
system: Application to hydrogen-bond breaking in water. J. Chem. 
Phys., 109:1125-1133, 1998. 

[34] P. G. Bolhuis, C. Dellago, and D. Chandler. Sampling ensembles of 
deterministic transition pathways. Faraday Discuss., 110:421-436, 1998. 

[35] C. Dellago, P. G. Bolhuis, and D. Chandler. On the calculation of 
reaction rate constants in the transition path ensemble. J. Chem. Phys., 
110:6617-6625, 1998. 

[36] D. Chandler. Barrier crossings: classical theory of rare but important 
events. In B.J. Berne, G. Ciccotti, and D.F. Coker, editors, Classical 
and Quantum Dynamics in Condensed Phase Simulations, pages 3-23. 
World Scientific Press, Singapore, 1998. 

[37] R. Olender and R. Elber. Calculation of classical trajectories with a very 
large time step: Formalism and numerical examples. J. Chem. Phys., 
105:9299-9315, 1996. 

[38] R. Elber, J. Meller, and R. Olender. Stochastic path approach to com- 
pute atomically detailed trajectories: application to the folding of C 
peptide. J. Phys. Chem. B, 103:899-911, 1999. 

[39] A. Ghosh, R. Elber, and H. A. Scheraga. An atomically detailed study of 
the folding pathways of protein a with the stochastic difference equation. 
Proc. Nat. Acad. Sci., 99:10394-10398, 2002. 

[40] T. B. Woolf. Path corrected functionals of stochastic trajectories: to- 
wards relative free energy and reaction coordinate calculations. Chem. 
Phys. Lett, 289:433-441, 1998. 

[41] D. M. Zuckerman and T. B. Woolf. Dynamic reaction paths and rates 
through importance-sampled stochastic dynamics. J. Chem. Phys., 
111:9475-9484, 1999. 

[42] D. M. Zuckerman and T. B. Woolf. Efficient dynamic importance sam- 
pling of rare events in one dimension. Phys. Rev. E, 63:016702, 2001. 

[43] D. M. Zuckerman and T. B. Woolf. Rapid determination of multiple 
reaction pathways in molecular systems: The soft ratcheting algorithm. 
www.arXiv.org: |physics/0209098| 2002. 

[44] O. Mazonka, C. Jarzyhski, and J. Blocki. Computing probabilities of 
very rare events for Langevin processes: a new method based on impor- 
tance sampling. Nuc. Phys. A, 641:335-354, 1998. Erratum: Nuc. Phys. 
A, 650, 499-500 (1999). 

[45] P. Eastman, N. Gr0nbech- Jensen, and S. Doniach. Simulation of protein 
folding by reaction path annealing, (unpublished), 2000. 



[46; 



G. A. Huber and S. Kim. Weighted-ensemble Brownian dynamics sim- 
ulations for protein association reactions. Biophys. J., 70:97-110, 1996. 



[58 



A. L. Borovinskiy and A. Y. Grosberg. Design of toy pro- 
teins capable to rearrange conformations in a mechanical fashion, 
www. ar Xiv. org ]cond-mat /02 1 2 1 24| (2002). 

A. Z. Panagiotopoulos and S. K. Kumar. Large lattice discretization 
effects on the phase coexistence of ionic fluids. Phys. Rev. Lett., 83:2981- 
2984, 1999. 

A. Z. Panagiotopoulos. On the equivalence of continuum and lattice 
models for fluids. J. Chem. Phys., 112:7132-7137, 2000. 

Y. Ueda, H. Taketomi, and N. Go. Studies on protein folding, unfolding 
and fluctuations by computer simulation. II. a three-dimensional lattice 
model of lysozyme. Biopolymers, 17:1531-1548, 1978. 

J. Skolnick, A. Kolinski, and R. Yaris. Monte Carlo simulations of the 
folding of beta-barrel globular proteins. Proc. Nat. Acad. Sci., 85:5057- 
5061, 1988. 

M. Vieth, A. Kolinski, C. L. Brooks, and J. Skolnick. Prediction of the 
folding pathways and structure of the GCN4 leucine zipper. J. Molec. 
Bio., 237:361-367, 1994. 

A. Kolinski and J. Skolnick. Lattice models of protein folding, dynamics, 
and thermodynamics. Chapman and Hall, New York, 1996. 

K. F. Lau and K. A. Dill. A lattice statistical mechan-ics model of the 
conformation and sequence spaces of proteins. Macromolecules, 22:3986- 
3997, 1989. 

H. S. Chan and K. A. Dill. Intrachain loops in polymers: effects of 
excluded volume. J. Chem. Phys., 90:492-509, 1989. 

K. F. Lau and K. A. Dill. Theory for protein mutability and biogenesis. 
Proc. Nat. Acad. Set., 87:638-642, 1990. 

D. G. Covell and R. L. Jernigan. Conformations of folded proteins in 
restricted spaces. Biochem., 29:3287-3294, 1990. 

Y. Ueda, H. Taketomi, and N. Go. Studies on protein folding, unfolding 
and fluctuations by computer simulation, i. the effects of specie amino 
acid sequences represented by specie inter-unit interactions. Int. J. Pep- 
tide Protein Res., 7:445-459, 1975. 



[59] C. Clementi, H. Nymeyer, and J. N. Onuchic. Topological and ener- 
getic factors: what determines the structural details of the transition 
state ensemble and "en route" intermediates for protein folding? an 
investigation for small globular proteins. J. Molec. Bio., 298:937-953, 
2000. 



[60] C. Clementi, P. A. Jennings, and J. N. Onuchic. How native-state topol- 
ogy affects the folding of dihydrofolate reductase and interleukin-lbeta. 
P. Nat. Acad. Sci, 97:5871-5876, 2000. 

[61] J. Shimada and E. I. Shakhnovich. The ensemble folding kinetics of 
protein G from an all- atom Monte Carlo simulation. P. Nat. Acad. Sci., 
99:11175-11180, 2002. 

[62] M. Levitt and A. Warshel. Computer simulation of protein folding. 
Nature, 253:694-698, 1975. 

[63] M. Levitt. A simplified representation of protein conformations for rapid 
simulation of protein folding. J. Molec. Bio., 104:59-107, 1976. 

[64] S. Miyazawa and R. L. Jernigan. Estimation of effective interresidue 
contact energies from protein crystal structures: quasi-chemical approx- 
imation. Macromol, 18:534-552, 1985. 

[65] S. Miyazawa and R. L. Jernigan. Residue-residue potentials with a 
favorable contact pair term and an unfavorable high packing density 
term, for simulation and threading. J. Molec. Bio., 256:623-644, 1996. 

[66] C. Wilson and S. Doniach. A computer model to dynamically simulate 
protein folding. Proteins, 6:193-209, 1989. 

[67] M. S. Friedrichs and P. G. Wolynes. Toward protein tertiary struc- 
ture recognition by means of associative memory hamiltonians. Science, 
246:371-373, 1989. 

[68] M. S. Friedrichs, R. A. Goldstein, and P. G. Wolynes. Generalized 
protein tertiary structure recognition using associative memory hamil- 
tonians. J. Molec. Bio., 222:1013-1034, 1991. 

[69] J. D. Honeycutt and D. Thirumalai. Metastability of the folded states 
of globular proteins. Proc. Nat. Acad. Set., 87:3526-3529, 1990. 

[70] Z. Guo, D. Thirumalai, and J. D. Honeycutt. Folding kinetics of pro- 
teins: a model study. J. Chem. Phys., 97:525-535, 1992. 

[71] J. D. Honeycutt and D. Thirumalai. The nature of folded states of 
globular proteins. Biopolymers, 32:695-709, 1992. 

[72] Y. Zhou, C. K. Hall, and M. Karplus. First-order disorder-to-order 
transition in an isolated homopolymer model. Phys. Rev. Lett., 77:2822- 
2825, 1996. 

[73] Y. Zhou and M. Karplus. Folding thermodynamics of a model three- 
helix bundle protein. P. Nat. Acad. Sci, 94:14429-14432, 1997. 

[74] Y. Zhou and M. Karplus. Folding of a model three-helix bundle protein: 
a thermodynamic and kinetic analysis. J. Molec. Bio., 293:917-951, 
1999. 



[75] A. V. Smith and C. K. Hall. Alpha-helix formation: discontinuous 
molecular dynamics on an intermediate-resolution protein model. Pro- 
teins, 44:344-360, 2001. 

[76] A. V. Smith and C. K. Hall. Assembly of a tetrameric alpha-helical bun- 
dle: computer simulations on an intermediate-resolution protein model. 
Proteins, 44:376-391, 2001. 

[77] C. Micheletti, J. R. Banavar, and A. Maritan. Confromations of proteins 
in equilibrium. Phys. Rev. Lett, 87:088102, 2001. 

[78] S. Tanaka and H. A. Scheraga. Medium- and long-ranged interaction 
parameters between amino acids fro predicting three-dimensional struc- 
tures of proteins. Macromolecules, 9:945-950, 1976. 

[79] A. Liwo, S. Oldziej, M. R. Pincus, R. J. Wawak, S. Rackovsky, and 
H. A. Scheraga. A united-residue force field for off-lattice protein- 
structure simulations. I. functional forms and parameters of long-range 
side-chain interaction potentials from protein crystal data. J. Comp. 
Chem., 18:849-873, 1997. 

[80] A. Liwo, M. R. Pincus, R. J. Wawak, S. Rackovsky, S. Oldziej, and H. A. 
Scheraga. A united-residue force field for off-lattice protein-structure 
simulations. II. parameterization of short-range interactions and deter- 
mination of weights of energy terms by z-score optimization. J. Comp. 
Chem., 18:874-887, 1997. 

[81] A. Liwo, C. Czaplewski, J. Pillardy, and H. A. Scheraga. Cumulant- 
based expressions for the multibody terms for the correlation between 
local and electrostatic interactions in the united-residue force field. J. 
Chem. Phys., 115:2323-2347, 2001. 

[82] A. Monge, R. A. Friesner, and B. Honig. An algorithm to generate 
low-resolution protein tertiary structures from knowledge of secondary 
structure. Proc. Nat. Acad. Set., 91:5027-5029, 1994. 

[83] J. R. Gunn, A. Monge, R. A. Friesner, and C. H. Marshall. Hierarchical 
algorithm for computer modeling of protein tertiary structure: folding 
of myoglobin to 6.2-angstrom resolution. J. Phys. Chem., 98:702-711, 
1994. 

[84] A. Monge, E. J. P. Lathrop, J. R. Gunn, P. S. Shenkin, and R. A. 
Friesner. Computer modeling of protein folding: conformational and en- 
ergetic analysis of reduced and detailed protein models. J. Molec. Bio., 
247:995-1012, 1995. 

[85] D. Tobi, G. Shafran, N. Linial, and R. Elber. On the design and analysis 
of protein folding potentials. Proteins, 40:71-85, 2000. 

[86] R. Srinivasan and G.D. Rose. LINUS: A hierarchic procedure to predict 
the fold of a protein. Proteins, 22:81-99, 1995. 



[87] I. Bahar, B. Erman, T. Halilglu, and R. L. Jernigan. Efficient charac- 
terization of collective motions and interresidue correlations in proteins 
by low-resolution simulations. Biochem., 36:13512-13523, 1997. 

[88] R. L. Jernigan and I. Bahar. Structure-derived potentials and protein 
simulations. Curr. Op. Struct. Bio., 6:195-209, 1996. 

[89] M. M. Tirion. Large amplitude elastic motions in proteins from a single- 
parameter, atomic analysis. Phys. Rev. Lett., 77:1905-1908, 1996. 

[90] I. Bahar, A. R. Atilgan, and B. Erman. Direct evaluation of thermal 
fluctuations in proteins using a single-parameter harmonic model. Fold. 
Design, 2:173-181, 1997. 

[91] T. Haliloglu, I. Bahar, and B. Erman. Gaussian dynamics of folded 
proteins. Phys. Rev. Lett, 79:3090-3093, 1997. 

[92] I. Bahar, A. R. Atilgan, M. C. Demirel, and B. Erman. Vibrational 
dynamics of folded proteins: significance of slow and fast motions in 
relation to function and stability. Phys. Rev. Lett., 80:2733-2736, 1998. 

[93] I. Bahar. Dynanics of proteins and biomolecular complexes: inferring 
functional motions from structure. Rev. Chem. Eng., 15:319-347, 1999. 

[94] K. Hinsen. Analysis of domain motions by approximate normal mode 
calculations. Proteins, 33:417-429, 1998. 

[95] M. P. Allen and D. J. Tildesley. Computer Simulation of Liquids. Oxford 
University Press, Oxford, 1987. 

[96] D. A. McQuarrie. Statistical Mechanics. Harper and Row, New York, 
New York, 1976. 

[97] M. H. Zaman, M. Shen, R. S. Berry, and K. F. Freed. Computer sim- 
ulation of met-enkephalin using explicit atom and united atom poten- 
tials: similarities, differences, and suggestions for improvement. J. Phys. 
Chem. B, 107:1685-1691, 2003. 



